set more off

* the data used in this paper may not be publicly released and aspects of the data cleaning files may violate terms of our data use agreeements.
* our code takes a final analysis dataset as given, and provides code of all actual analyses in the paper.
* this file reflects analysis for period 1 (2002-2005), bottled water.

*************************************************************************************************************************************
* ANALYSIS CONSTRUCTION 2002-2005
*************************************************************************************************************************************

	**************************************************************************************************
	* THIS IS WHERE LEVEL IS CHOSEN
	**************************************************************************************************
	
	* level is county (0) or store (1)
	local level=1

if `level'==1 {
	use period1_store_treatments, clear
	tsset, clear
	merge 1:1 storeid statadate using period1_store_sales
}
else {
	use period1_county_treatments, clear
	tsset, clear
	merge 1:1 countyid statadate using period1_county_sales
}

keep if _merge==3
	* merge=2 is months with no possible treatments (i.e. dec 2005)
	* merge=1 is stores or time periods without sales data

drop _merge

destring countyid, gen(countyid_num)
if `level'==1 {
	encode storeid, gen(storeid_num)
}

sort statadate
egen week=group(statadate)

if `level'==1 {
	tsset storeid_num week
}
else {
	tsset countyid_num week
}

******************************************************************************************************
* EVENT STUDY
*****************************************************************************************************

* note that these are point-wise confidence intervals

if `level'==1 { 
	tsset storeid_num week


foreach z in "water" {

preserve

sum expend_`z'
gen mean=r(mean)
replace expend_`z'=expend_`z'-mean
replace expend_`z'=expend_`z'/1000

local win = 10
forvalues x=0/`win' {
	gen expend_`z'_l`x'=l`x'.expend_`z'
}
forvalues x=1/`win' {
	gen expend_`z'_f`x'=f`x'.expend_`z'
}

*keep if hlandfall_false==1
keep if hlandfall==1

keep expend_`z'_*

forvalues x=0/`win' {
	gen expend2_`z'_l`x'=expend_`z'_l`x'
}
forvalues x=1/`win' {
	gen expend2_`z'_f`x'=expend_`z'_f`x'
}

collapse (sem) expend2_`z'_* (mean) expend_`z'_*

duplicates drop
forvalues x=1/`win' {
	local y=`x'+50
	rename expend_`z'_f`x' expend_`z'_`y'
	rename expend2_`z'_f`x' expend2_`z'_`y'
	
}
forvalues x=0/`win' {
	local y=50-`x'
	rename expend_`z'_l`x' expend_`z'_`y'
	rename expend2_`z'_l`x' expend2_`z'_`y'
}

gen id=1
reshape long expend_`z'_ expend2_`z'_, i(id) j(time)
drop id
replace time=time-50

gen ci=expend2_`z'_*1.96
gen ul=expend_`z'+ci
gen ll=expend_`z'-ci

twoway (rarea ul ll time, color(gs14)) (connected expend_`z'_ time, msize(small)) (pcarrowi 0 0 4 0, msize(zero)), scheme(s2mono) plotr(m(zero)) ytitle("") ylabel(0(1)4) xlabel(-10(2)10) xtitle("Weeks Before and After Storm Landfall") legend(off) title("Water Sales - 2002-2005") saving(Fig2A_false, replace)

restore

}

}

****************************************************************************************************
* ADDITIONAL PROCESSING FOR REGRESSION ANALYSIS
****************************************************************************************************

* merge in weather at store level 
* some winter weeks may be missing weather data. we do not use winter weeks, so do not worry about this
* merge==2 is weeks with no sales 
if `level' == 1 {
	merge 1:1 storeid statadate using period1_store_weather
}
else {
	merge 1:1 countyid statadate using period1_county_weather
}
drop if _merge==2
drop _merge

*create year dummies
tab year, gen(years)

* generate logged depvars - NOTE ADD $1 and log
foreach code in "water" "soda" "dietsoda" "milk" "fruitdrinks" "apple" "orange" {
	gen lexpend_`code'=log(expend_`code'+1)
	}

* create month and day variables
capture drop month
capture drop day
gen month=month(statadate)
gen day=day(statadate)

* define first week of hurricane season 
if `level' == 1 {

gen first=((year==2002&month==6&day==1)|(year==2003&month==5&day==31)|(year==2004&month==5&day==29)|(year==2005&month==5&day==28))

* define week of first hurricane of that year's hurricane season. 
* note set to missing if that store is impacted by a hurricane sometime during that year.
bys year: egen temp=min(statadate) if hlandfall==1
bys year: egen temp2 =min(temp)
format temp2 %td
gen firsth=(temp2==statadate)
drop temp temp2
bys storeid year: egen temp1=max(htreat_before)
bys storeid year: egen temp2=max(htreat_after)
bys storeid year: egen temp3=max(hlandfall)
egen ever=rowmax(temp1 temp2 temp3)
replace firsth=. if ever==1
drop temp1 temp2 temp3

}

* several demos and risk factors already in the store level data. we add factors that are not included.
if `level' == 0 {
	merge m:1 countyid using county_demos
	keep if _m==3
	drop _merge
}
	
* merge in risk factors at county level
merge m:1 countyid using county_coastal
drop if _merge==2
gen coastal_county=(coastal=="C")
drop coastal _merge

if `level' == 1 {
	gen coastal_store=(distance<=5)
}

if `level' == 1 {
	merge m:1 countyid statadate using period1_county_fema
}
else {
	merge 1:1 countyid statadate using period1_county_fema
}	
drop if _merge==2 
gen fema_county_week=(fema==1)
drop fema _merge

if `level' == 0 {
	rename pctblack1 pblack
	rename medianhhinc hinc00
	rename pctbachelorsormore pcol
}

foreach var in pblack hinc00 pcol {
	quietly sum `var', detail
	scalar m`var'=r(p50)
	gen hi`var'=(`var'>m`var') if `var'!=.
	gen b_inter_hi`var'=hi`var'*htreat_before
	gen d_inter_hi`var'=hi`var'*hlandfall
	gen a_inter_hi`var'=hi`var'*htreat_after
}

if `level' == 1 {

	foreach var in pblack hinc00 pcol {
		gen first_inter_hi`var'=hi`var'*first
		gen firsth_inter_hi`var'=hi`var'*firsth
	}
}

foreach rvar in coastal_county fema_county_week {
	gen b_inter_`rvar'=`rvar'*htreat_before
	gen d_inter_`rvar'=`rvar'*hlandfall
	gen a_inter_`rvar'=`rvar'*htreat_after
	gen b_inter_`rvar'_evac=`rvar'*htreat_before_evac
	gen d_inter_`rvar'_evac=`rvar'*hlandfall_evac
	gen a_inter_`rvar'_evac=`rvar'*htreat_after_evac
	gen b_inter_`rvar'_nonevac=`rvar'*htreat_before_nonevac
	gen d_inter_`rvar'_nonevac=`rvar'*hlandfall_nonevac
	gen a_inter_`rvar'_nonevac=`rvar'*htreat_after_nonevac
}

if `level' == 1 {

	foreach rvar in coastal_county fema_county_week {

		gen first_inter_`rvar'=`rvar'*first
		gen firsth_inter_`rvar'=`rvar'*firsth
	}
}

if `level' == 1 {

	foreach r2var in coastal_store {
		gen b_inter_`r2var'=`r2var'*htreat_before
		gen d_inter_`r2var'=`r2var'*hlandfall
		gen a_inter_`r2var'=`r2var'*htreat_after
		gen b_inter_`r2var'_evac=`r2var'*htreat_before_evac
		gen d_inter_`r2var'_evac=`r2var'*hlandfall_evac
		gen a_inter_`r2var'_evac=`r2var'*htreat_after_evac
		gen b_inter_`r2var'_nonevac=`r2var'*htreat_before_nonevac
		gen d_inter_`r2var'_nonevac=`r2var'*hlandfall_nonevac
		gen a_inter_`r2var'_nonevac=`r2var'*htreat_after_nonevac
		gen first_inter_`r2var'=`r2var'*first
		gen firsth_inter_`r2var'=`r2var'*firsth
	}
}


* merge in sales by size
if `level' == 1 {

	merge 1:1 storeid statadate using period1_store_revbysize
	drop if _merge==2
	* merge 1 is master data with missing sales weeks for squaring
	* merge 2 is revbysize data from stores not in final sample 
	drop _merge

	foreach num of numlist 16906 16912 128001 {
		replace expend_`num'=0 if expend_`num'==. & expend_water!=.
		gen lexpend_`num'=log(expend_`num'+1)
	}
}

* assign fixed effects level
if `level'==1 {
	local FE "storeid_num"	 
}
else {
	local FE "countyid_num"
}

* generate impact leads
if `level'==1 {
	tsset storeid_num week 
}
else {
	tsset countyid_num week
}
gen hlandfall_f1=f1.hlandfall
gen hlandfall_f2=f2.hlandfall

***********************************************************************************************************************************
* REGRESSIONS
***********************************************************************************************************************************

* early season stockpiling 

if `level'==1 {

preserve

	* this analysis must be done before full restrictions to hurricane season only.
	* here, we restrict analysis to april - november to allow some pre-season data but to otherwise keep analysis comparable.

keep if inrange(month,5,11)
tab month, gen(months)

xtreg lexpend_water first prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water first first_inter_hipblack prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water first first_inter_hihinc00 prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water first first_inter_hipcol prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water first first_inter_coastal_county prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)

xtreg lexpend_water firsth prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water firsth firsth_inter_hipblack prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water firsth firsth_inter_hihinc00 prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water firsth firsth_inter_hipcol prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)
xtreg lexpend_water firsth firsth_inter_coastal_county prcp tmean tmax tmin months1-months6 years1-years3,  fe i(`FE') vce(cluster countyid_num)

restore

}

* MAIN ANALYSIS

* keep just hurricane season.

keep if inrange(month,6,11)
tab month, gen(months)


if `level'==1 {

* period 1 main regressions + other beverages

foreach code in "water" "soda" "dietsoda" "apple" "orange" {
	xtreg lexpend_`code' htreat_before hlandfall htreat_after prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)
}

* period1 regressions with interactions. with forecast, weather, simplified treatment definitions.

xtreg lexpend_water htreat_before b_inter_hipblack hlandfall d_inter_hipblack htreat_after a_inter_hipblack prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)

xtreg lexpend_water htreat_before b_inter_hihinc00 hlandfall d_inter_hihinc00 htreat_after a_inter_hihinc00 prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)

xtreg lexpend_water htreat_before b_inter_hipcol hlandfall d_inter_hipcol htreat_after a_inter_hipcol prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)

xtreg lexpend_water htreat_before b_inter_coastal_county hlandfall d_inter_coastal_county htreat_after a_inter_coastal_county prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)

xtreg lexpend_water htreat_before b_inter_fema_county_week hlandfall d_inter_fema_county_week htreat_after a_inter_fema_county_week fema_county_week prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)


* period1 analysis by common sizes

foreach num of numlist 16906 16912 128001 {
	xtreg lexpend_`num' htreat_before hlandfall htreat_after prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)
}

* period1 analysis by evacuation status

foreach status in "evac" "nonevac" {
	xtreg lexpend_water htreat_before_`status' hlandfall_`status' htreat_after_`status' prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)
	xtreg lexpend_water htreat_before_`status' b_inter_coastal_county_`status' hlandfall_`status' d_inter_coastal_county_`status' htreat_after_`status' a_inter_coastal_county_`status' prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)
}

* threatened near miss vs. threatened followed by strike. 

	* is this store week impacted by hurricane this week, next week, or week after that?
	egen temp2=rowmax(hlandfall hlandfall_f1 hlandfall_f2)
	* is this store week threatened by a hurricane this week, but not impacted by hurricane this week, next week, or week after that?
	gen nearmiss = (temp2==0 & htreat_before==1)
	* is this store week threatened by a hurricane this week, and subsequently impacted by hurricane this week, next week, or week after that?
	gen threatened_hit = (temp2==1 & htreat_before==1)

	* how many store weeks are hit by a storm but not threatened?
	bys storeid_num year: egen temp3=max(threatened_hit)
	count if hlandfall==1 & temp3==0
	
	* summary stats
	tab htreat_before
	tab nearmiss
	tab threatened_hit

	* generate impact leads
	tsset storeid_num week 
	gen nearmiss_l1=l1.nearmiss
	gen nearmiss_l2=l2.nearmiss

	bys week: egen temp33=max(hlandfall)
	bys week: egen temp44=max(htreat_after)

	* creating simulated impacts
	gen nearmiss_during = (nearmiss==1 | nearmiss_l1==1 | nearmiss_l2==1) & (temp33==1)
	gen nearmiss_after = (nearmiss==1 | nearmiss_l1==1 | nearmiss_l2==1) & (temp44==1)

	xtreg lexpend_water htreat_before hlandfall htreat_after prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)
	xtreg lexpend_water nearmiss nearmiss_during nearmiss_after prcp tmean tmax tmin months1-months5 years1-years3 ,  fe i(`FE') vce(cluster countyid_num)
	xtreg lexpend_water threatened_hit hlandfall htreat_after prcp tmean tmax tmin months1-months5 years1-years3 ,  fe i(`FE') vce(cluster countyid_num)

}

else {

* county-level main regression
xtreg lexpend_water htreat_before hlandfall htreat_after prcp tmean tmax tmin months1-months5 years1-years3,  fe i(`FE') vce(cluster countyid_num)

}

STOP
STOP
STOP
STOP



